Hemos elegido el dataset “Wine quality dataset”, un dataset de kaggle cuya fuente es el repositorio UCI Machine Learning (https://www.kaggle.com/uciml/red-wine-quality-cortez-et-al-2009; https://archive.ics.uci.edu/ml/datasets/Wine+Quality). Se trata de un dataset de 6497 observaciones, 11 atributos y 1 marca de clase.
Resumen de datos: más precisamente son 2 datasets separados, para tipos tinto y blanco del vinho verde, el de vinos tintos contiene 1599 observaciones, y el de vinos blancos - 4898 observaciones, con que tenemos opción de trabajar con datasets separadamente o bien unir todas las observaciones y obtener la variebilidad natural tinto/blanco.
Las clases de dataset corresopnden a la calidad del vino basandose en sus características físicas y químicas (sensory data), en una escala de 0 a 10. El paradigma de calidad y las características son facilmente interpretables y el ámbito del tema es generalmente conocido pues tiene un alto potencial de dar un resultado divulgativo.
En cuanto a la finalidad del estudio, sería interesante jugar con los datos para estimar que caracteristicas influyen a la calidad y/o representan el color. ¿Las caracteristicas de calidad varian segun el color del vino? Como podemos distribuir la calidad para entenderlo empiricamente -qué agrupaciones podemos obtener y cual es el número de “categorías de calidad” más óptimo según los algoritmos?, o bien si podemos estudiar los vinos partiendo de otras variables.
# Cargamos las librerías
library(ggplot2) # visualization
library(patchwork)
library(tidyverse)
library(corrplot)
library(factoextra) # clustering
library(Hmisc) # impute
library(arules) # discretize
library(DMwR) # smote
library(DescTools) # box-cox
Cargamos los dos datasets y concatenamos los datos de vinos tintos y blancos:
# Carga del dataset
red <- read.csv('https://archive.ics.uci.edu/ml/machine-learning-databases/wine-quality/winequality-red.csv', header = TRUE, sep = ";")
white <- read.csv('https://archive.ics.uci.edu/ml/machine-learning-databases/wine-quality/winequality-white.csv', header = TRUE, sep = ";")
# Nombres de los atributos
names(red) <- c("fixed acidity", "volatile acidity", "citric acid", "residual sugar", "chlorides", "free sulfur dioxide", "total sulfur dioxide" , "density", "pH", "sulphates", "alcohol", "quality")
names(white) <- c("fixed acidity", "volatile acidity", "citric acid", "residual sugar", "chlorides", "free sulfur dioxide", "total sulfur dioxide" , "density", "pH", "sulphates", "alcohol", "quality")
# Creamos las columnas de color
red['color'] <- 'red'
white['color'] <- 'white'
# Fusionamos los dos datasets
vinos <- rbind(white,red)
# Factorizamos la variable color
vinos$color <- as.factor(vinos$color)
# Verificamos la estructura del conjunto de datos
str(vinos)
## 'data.frame': 6497 obs. of 13 variables:
## $ fixed acidity : num 7 6.3 8.1 7.2 7.2 8.1 6.2 7 6.3 8.1 ...
## $ volatile acidity : num 0.27 0.3 0.28 0.23 0.23 0.28 0.32 0.27 0.3 0.22 ...
## $ citric acid : num 0.36 0.34 0.4 0.32 0.32 0.4 0.16 0.36 0.34 0.43 ...
## $ residual sugar : num 20.7 1.6 6.9 8.5 8.5 6.9 7 20.7 1.6 1.5 ...
## $ chlorides : num 0.045 0.049 0.05 0.058 0.058 0.05 0.045 0.045 0.049 0.044 ...
## $ free sulfur dioxide : num 45 14 30 47 47 30 30 45 14 28 ...
## $ total sulfur dioxide: num 170 132 97 186 186 97 136 170 132 129 ...
## $ density : num 1.001 0.994 0.995 0.996 0.996 ...
## $ pH : num 3 3.3 3.26 3.19 3.19 3.26 3.18 3 3.3 3.22 ...
## $ sulphates : num 0.45 0.49 0.44 0.4 0.4 0.44 0.47 0.45 0.49 0.45 ...
## $ alcohol : num 8.8 9.5 10.1 9.9 9.9 10.1 9.6 8.8 9.5 11 ...
## $ quality : int 6 6 6 6 6 6 6 6 6 6 ...
## $ color : Factor w/ 2 levels "red","white": 2 2 2 2 2 2 2 2 2 2 ...
Se observa que tenemos todas las variables numéricas menos el atributo de color.
Descripción de variables: son las características físicas o químicas del vino más su calidad y su color/tipo (tinto o blanco)
fixed acidity : acidez fija (g/l), v. continua
volatile acidity : acidez volatil (g/l), v. continua
citric acid : acido citrico (g/l), v. continua
residual sugar : azucar residual (g/l), v. continua
chlorides : cloruros (g/l), v. continua
free sulfur dioxide : dioxido de azufre libre (mg/l), v. continua
total sulfur dioxide: dioxido de azufre total (mg/l), v. continua
density : densidad (g/l), v. continua
pH : pH, v. continua
sulphates : sulfatos (g/l), v. continua
alcohol : concentracción de alcohol (%%), v. continua
quality : calidad, v. discreta
color : color, v. cualitativa dicotómica
Veamos como son algunas de las observaciones:
rbind(head(vinos,5), tail(vinos,5))
## fixed acidity volatile acidity citric acid residual sugar chlorides
## 1 7.0 0.270 0.36 20.7 0.045
## 2 6.3 0.300 0.34 1.6 0.049
## 3 8.1 0.280 0.40 6.9 0.050
## 4 7.2 0.230 0.32 8.5 0.058
## 5 7.2 0.230 0.32 8.5 0.058
## 6493 6.2 0.600 0.08 2.0 0.090
## 6494 5.9 0.550 0.10 2.2 0.062
## 6495 6.3 0.510 0.13 2.3 0.076
## 6496 5.9 0.645 0.12 2.0 0.075
## 6497 6.0 0.310 0.47 3.6 0.067
## free sulfur dioxide total sulfur dioxide density pH sulphates
## 1 45 170 1.00100 3.00 0.45
## 2 14 132 0.99400 3.30 0.49
## 3 30 97 0.99510 3.26 0.44
## 4 47 186 0.99560 3.19 0.40
## 5 47 186 0.99560 3.19 0.40
## 6493 32 44 0.99490 3.45 0.58
## 6494 39 51 0.99512 3.52 0.76
## 6495 29 40 0.99574 3.42 0.75
## 6496 32 44 0.99547 3.57 0.71
## 6497 18 42 0.99549 3.39 0.66
## alcohol quality color
## 1 8.8 6 white
## 2 9.5 6 white
## 3 10.1 6 white
## 4 9.9 6 white
## 5 9.9 6 white
## 6493 10.5 5 red
## 6494 11.2 6 red
## 6495 11.0 6 red
## 6496 10.2 5 red
## 6497 11.0 6 red
Veamos las estadísticas básicas:
summary(vinos)
## fixed acidity volatile acidity citric acid residual sugar
## Min. : 3.800 Min. :0.0800 Min. :0.0000 Min. : 0.600
## 1st Qu.: 6.400 1st Qu.:0.2300 1st Qu.:0.2500 1st Qu.: 1.800
## Median : 7.000 Median :0.2900 Median :0.3100 Median : 3.000
## Mean : 7.215 Mean :0.3397 Mean :0.3186 Mean : 5.443
## 3rd Qu.: 7.700 3rd Qu.:0.4000 3rd Qu.:0.3900 3rd Qu.: 8.100
## Max. :15.900 Max. :1.5800 Max. :1.6600 Max. :65.800
## chlorides free sulfur dioxide total sulfur dioxide
## Min. :0.00900 Min. : 1.00 Min. : 6.0
## 1st Qu.:0.03800 1st Qu.: 17.00 1st Qu.: 77.0
## Median :0.04700 Median : 29.00 Median :118.0
## Mean :0.05603 Mean : 30.53 Mean :115.7
## 3rd Qu.:0.06500 3rd Qu.: 41.00 3rd Qu.:156.0
## Max. :0.61100 Max. :289.00 Max. :440.0
## density pH sulphates alcohol
## Min. :0.9871 Min. :2.720 Min. :0.2200 Min. : 8.00
## 1st Qu.:0.9923 1st Qu.:3.110 1st Qu.:0.4300 1st Qu.: 9.50
## Median :0.9949 Median :3.210 Median :0.5100 Median :10.30
## Mean :0.9947 Mean :3.219 Mean :0.5313 Mean :10.49
## 3rd Qu.:0.9970 3rd Qu.:3.320 3rd Qu.:0.6000 3rd Qu.:11.30
## Max. :1.0390 Max. :4.010 Max. :2.0000 Max. :14.90
## quality color
## Min. :3.000 red :1599
## 1st Qu.:5.000 white:4898
## Median :6.000
## Mean :5.818
## 3rd Qu.:6.000
## Max. :9.000
Todas los columnas parecen ser bastante limpias, no obstante aquí se observa la heterogeneidad de las variables (por ejemplo, cloruros que no superan 0.611 g/l y dioxido de sulfuro total que “alcanza” 440 mg/l: dado que tienen las unidades distintas se produce esta brecha).
En cuanto a la calidad, el valor mínimo es 3 y el máximo es 9. Por lo que 1-10 es una escala “común”, y la escala real sería 3-9. Comprobamos la distribución de calidad:
table(vinos$quality)
##
## 3 4 5 6 7 8 9
## 30 216 2138 2836 1079 193 5
Con el atributo “quality” podemos tener 7 marcas de clase diferentes, aunque puede ser conveniente agruparlo en una variable discreta como vemos más adelante.
Distribución de valores únicos de atributos:
apply(vinos,2, function(x) length(unique(x)))
## fixed acidity volatile acidity citric acid
## 106 187 89
## residual sugar chlorides free sulfur dioxide
## 316 214 135
## total sulfur dioxide density pH
## 276 998 108
## sulphates alcohol quality
## 111 111 7
## color
## 2
Algunos atributos tienen la distribución bastante dispersa, que, si fuera necesario, podríamos agrupar en intervalos.
Comprobamos si hay valores nulos en en dataset
any(is.na(vinos))
## [1] FALSE
any(vinos=="")
## [1] FALSE
A partir de las estadísticas se ve que todas las variables tienen un mínimo distinto del cero menos la variable “citric acid”, y podría tomar un 0 como un valor desconocido.
Veamos la distribución de “citric acid”:
table(vinos$`citric acid`)
##
## 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1 0.11 0.12 0.13 0.14
## 151 40 56 32 41 25 30 34 37 42 49 16 46 35 48
## 0.15 0.16 0.17 0.18 0.19 0.2 0.21 0.22 0.23 0.24 0.25 0.26 0.27 0.28 0.29
## 42 42 43 71 69 95 99 131 108 232 163 257 236 301 244
## 0.3 0.31 0.32 0.33 0.34 0.35 0.36 0.37 0.38 0.39 0.4 0.41 0.42 0.43 0.44
## 337 230 289 208 249 150 197 153 136 129 146 98 124 52 86
## 0.45 0.46 0.47 0.48 0.49 0.5 0.51 0.52 0.53 0.54 0.55 0.56 0.57 0.58 0.59
## 68 70 56 62 283 55 38 40 30 32 23 30 22 30 14
## 0.6 0.61 0.62 0.63 0.64 0.65 0.66 0.67 0.68 0.69 0.7 0.71 0.72 0.73 0.74
## 15 11 15 14 15 15 21 9 18 9 5 10 6 8 45
## 0.75 0.76 0.78 0.79 0.8 0.81 0.82 0.86 0.88 0.91 0.99 1 1.23 1.66
## 1 3 3 3 2 2 2 1 1 2 1 6 1 1
Observación “0” aparece 151 veces que supone apr. 2% de la distribución y es comparable con algunas otras frecuencias, por lo que puede ser valor real (=“no siempre se añade el acido citrico”) y no perdido o nulo.
Como hemos visto, hay variables que tienen distintas unidades de medidas (hablando de variables de misma naturaleza), podemos reducirlas a las mismas medidas (quedamos con g/l).
Hay dos variables que usan otras unidades: “free/total sulfur dioxide”, que en g/l tendrían la distribución parecida a la de “chlorides”.
# Cambio de unidades de mg/l a g/l
vinos$`free sulfur dioxide` <- vinos$`free sulfur dioxide` * 0.001
vinos$`total sulfur dioxide` <- vinos$`total sulfur dioxide` * 0.001
Visualizamos de nuevo las estadísticas:
summary(vinos)
## fixed acidity volatile acidity citric acid residual sugar
## Min. : 3.800 Min. :0.0800 Min. :0.0000 Min. : 0.600
## 1st Qu.: 6.400 1st Qu.:0.2300 1st Qu.:0.2500 1st Qu.: 1.800
## Median : 7.000 Median :0.2900 Median :0.3100 Median : 3.000
## Mean : 7.215 Mean :0.3397 Mean :0.3186 Mean : 5.443
## 3rd Qu.: 7.700 3rd Qu.:0.4000 3rd Qu.:0.3900 3rd Qu.: 8.100
## Max. :15.900 Max. :1.5800 Max. :1.6600 Max. :65.800
## chlorides free sulfur dioxide total sulfur dioxide
## Min. :0.00900 Min. :0.00100 Min. :0.0060
## 1st Qu.:0.03800 1st Qu.:0.01700 1st Qu.:0.0770
## Median :0.04700 Median :0.02900 Median :0.1180
## Mean :0.05603 Mean :0.03053 Mean :0.1157
## 3rd Qu.:0.06500 3rd Qu.:0.04100 3rd Qu.:0.1560
## Max. :0.61100 Max. :0.28900 Max. :0.4400
## density pH sulphates alcohol
## Min. :0.9871 Min. :2.720 Min. :0.2200 Min. : 8.00
## 1st Qu.:0.9923 1st Qu.:3.110 1st Qu.:0.4300 1st Qu.: 9.50
## Median :0.9949 Median :3.210 Median :0.5100 Median :10.30
## Mean :0.9947 Mean :3.219 Mean :0.5313 Mean :10.49
## 3rd Qu.:0.9970 3rd Qu.:3.320 3rd Qu.:0.6000 3rd Qu.:11.30
## Max. :1.0390 Max. :4.010 Max. :2.0000 Max. :14.90
## quality color
## Min. :3.000 red :1599
## 1st Qu.:5.000 white:4898
## Median :6.000
## Mean :5.818
## 3rd Qu.:6.000
## Max. :9.000
El dataset parece tener outliers ya que en muchas variables la dfirenecia entre tercer quantil y el máximo es considerable.
Visualizamos los boxplots de todas las columnas:
bp1 <- vinos %>%
gather(Attributes, values, c(1:11)) %>%
ggplot(aes(x=Attributes, y=values, fill=Attributes)) + geom_boxplot(show.legend=FALSE) + coord_flip()
bp1
Se observa que la variable “residual sugar” tiene unos outliers muy distantes. Otros atributos tienen una escala diferente por lo que procedemos a visualizarlos sin “residual sugar” con distinta ampliación:
p1 <- vinos %>%
gather(Attributes, values, c(1:3,5:11)) %>%
ggplot(aes(x=Attributes, y=values, fill=Attributes)) + geom_boxplot(show.legend=FALSE) + coord_flip()
p2 <- vinos %>%
gather(Attributes, values, c(2:3,5:10)) %>%
ggplot(aes(x=Attributes, y=values, fill=Attributes)) + geom_boxplot(show.legend=FALSE) + coord_flip()
p3 <- vinos %>%
gather(Attributes, values, c(2:3,5:8,10)) %>%
ggplot(aes(x=Attributes, y=values, fill=Attributes)) + geom_boxplot(show.legend=FALSE) + coord_flip()
bp1 + p2 + p3 + plot_layout(nrow = 3)
Tenemos valores bastante anómalos, sobre todo en “citric acid”, “free sulfur dioxide”, “sulphates”, “volatile acidity”, “sulphates”, “chlorides”.
No obstante, son relativamente pocas observaciones por lo que se puede realizar una imputación por la mediana.
Por ello, remplazamos los valores extremos segun el estadístico de boxplot por NA:
for (x in c("residual sugar","citric acid", "free sulfur dioxide", "sulphates", "volatile acidity", "chlorides", "total sulfur dioxide", "fixed acidity", "density", "alcohol", "pH")) {
vinos[,x][vinos[,x] %in% (boxplot.stats(vinos[,x])$out) ] <- NA
}
Cantidad de outliers detectados:
sapply(vinos, function(vinos) sum(is.na(vinos)))
## fixed acidity volatile acidity citric acid
## 357 377 509
## residual sugar chlorides free sulfur dioxide
## 118 286 62
## total sulfur dioxide density pH
## 10 3 73
## sulphates alcohol quality
## 191 3 0
## color
## 0
Un ejemplo de observaciones con outliers:
vinos[is.na(vinos$alcohol),]
## fixed acidity volatile acidity citric acid residual sugar chlorides
## 3919 6.4 0.35 0.28 1.6 0.037
## 4504 5.8 0.61 NA 8.4 0.041
## 5551 NA 0.36 NA 7.5 0.096
## free sulfur dioxide total sulfur dioxide density pH sulphates
## 3919 0.031 0.113 0.98779 3.12 0.40
## 4504 0.031 0.104 0.99090 3.26 0.72
## 5551 0.022 0.071 0.99760 2.98 0.84
## alcohol quality color
## 3919 NA 7 white
## 4504 NA 7 white
## 5551 NA 5 red
Podemos ver que algunos valores atípicos se encuentran en las mismas observaciones.
Imputación por la mediana:
vinos[,c(1:11)] <- apply(vinos[,c(1:11)], 2, impute)
Visualizamos los atributos de nuevo comparando con los boxplots anteriores:
bp2 <- vinos %>%
gather(Attributes, values, c(1:11)) %>%
ggplot(aes(x=Attributes, y=values, fill=Attributes)) + geom_boxplot(show.legend=FALSE) + coord_flip()
bp1 + labs(subtitle = "Vinos original") + bp2 + labs(subtitle = "Vinos sin outliers") + plot_layout(nrow = 2)
O bien omitiendo residual sugar:
p2 <- vinos %>%
gather(Attributes, values, c(1:3,5:11)) %>%
ggplot(aes(x=Attributes, y=values, fill=Attributes)) + geom_boxplot(show.legend=FALSE) + coord_flip()
p1 + labs(subtitle = "Vinos original") + p2 + labs(subtitle = "Vinos sin outliers") + plot_layout(nrow = 2)
Aunque seguimos teniendo outliers segun estos boxplots en su definción más teórica, están más agrupados y, considerando que representan la variebilidad real del dataset, eliminamos los valores muy anómalos y demasiado dispersos que podían ser errores o inconsistencias. Por ello, hemos obtenido la distribución mucho menos sesgada.
Como se observa, la variable quality no está balanceada, y las clases que tienen pocas observaciones pueden presentar problemas en análisis así que es conveniente crear particiones con más observaciones. Por ello con el fin de equilibrar la marca de calidad y agruparlo de manera natural, se puede realizar la discretizacion de la variable.
Veamos de nuevo sus estadísticas:
summary(vinos$quality)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 3.000 5.000 6.000 5.818 6.000 9.000
barplot(table(vinos$quality), xlab="quality", main = "quality")
# Distribución de valores únicos
table(vinos$quality)
##
## 3 4 5 6 7 8 9
## 30 216 2138 2836 1079 193 5
Tenemos 7 marcas de calidad ordinales de 6471 observaciones. En principio, no sabemos en cuantas particiones podemos agrupar las observaciones.
Logicamente, por frecuencias no se puede crear 4 o más particiones equilibradas (6471 ÷ 4 = 1617, los números de observaciones de clases 5,6 son mayores). Por lo que podríamos crear 2 o 3 clases. Tenemos que se seguir la lógica real para agrupar la calidad, por lo que puede ser: alta/media/baja calidad o bien alta/no alta calidad. Visualizamos la distribución de clases segun el color:
ggplot(data = vinos,aes(x=color,fill=as.factor(quality)))+geom_bar(position="fill") + geom_hline(yintercept=0.333, alpha = 0.3) + geom_hline(yintercept=0.666, alpha = 0.3) + geom_hline(yintercept=0.5, alpha = 0.3, color='blue') + coord_flip()
Para poder trabajar posteriormente con una variable dicotómica de calidad, fijaremos el número de bins de 2, que representaria calidad alta/no alta.
Para poder determinar qué observaciones agrupamos en que clase de calidad, visualización de clases discretas segun si la particion se hace por igual frequencia / igual amplitud o clustering:
par(mfrow=c(2,2))
hist(vinos$quality, breaks = 20, main = "Equal Frequency")
abline(v = discretize(vinos$quality, breaks = 2, onlycuts = TRUE), col = "red")
hist(vinos$quality, breaks = 20, main = "Equal Width")
abline(v = discretize(vinos$quality, method = "interval", breaks = 2, onlycuts = TRUE), col = "red")
hist(vinos$quality, breaks = 20, main = "Clustering")
abline(v = discretize(vinos$quality, method = "cluster", breaks = 2, onlycuts = TRUE), col = "red")
Creamos las particiones como otro atributo para poder comparar las futuras agrupaciones con valores de quality. Aunque pueda ser redundante ahora, nos podría servir para ajustar las particiones de calidad en futuro.
vinos$`quality class`[vinos$quality<=5]="low"
vinos$`quality class`[vinos$quality>5]="high"
vinos$`quality class` <- as.factor(vinos$`quality class`)
# Creamos un factor ordinal
vinos$`quality class` <- factor(vinos$`quality class`, levels = c("low","high"))
# Estructura de la variable
table(vinos$`quality class`)
##
## low high
## 2384 4113
Finalmente, convertimos la calidad en factor también:
vinos$quality <- as.factor(vinos$quality)
# Variable factorizada
str(vinos$quality)
## Factor w/ 7 levels "3","4","5","6",..: 4 4 4 4 4 4 4 4 4 4 ...
Ahora tenemos un dataset de 6471 observaciones y 14 columnas, 3 de cuales se puede considerar clases / grupos de vinos reales: quality (7 niveles), quality class (2 niveles), color (2 niveles).
Si podemos considerar que los datos de clase de calidad están ahora balanceados, las observaciones de color no lo estan, por ello usamos el algortimo SMOTE para poder equilibrar las marcas de clase para tareas posteriores de análisis.
vinos <- SMOTE(color ~ ., vinos)
dim(vinos)
## [1] 11193 14
table(vinos$color)
##
## red white
## 4797 6396
table(vinos$quality)
##
## 3 4 5 6 7 8 9
## 56 373 4136 4760 1606 256 6
Con ello tenemos un dataset de 11193 observaciones con vinos tintos/blancos más balanceados gracias a las nuevas observaciones “sinteticas” y como se observa, la calidad tiene la distribución parecida a la anterior.
Volvemos a ver las estadicticas básicas de las variables:
summary(vinos)
## fixed acidity volatile acidity citric acid residual sugar
## Min. :4.500 Min. :0.0800 Min. :0.0400 Min. : 0.600
## 1st Qu.:6.599 1st Qu.:0.2500 1st Qu.:0.2502 1st Qu.: 1.900
## Median :6.999 Median :0.3000 Median :0.3100 Median : 2.490
## Mean :7.115 Mean :0.3381 Mean :0.3122 Mean : 4.471
## 3rd Qu.:7.624 3rd Qu.:0.4181 3rd Qu.:0.3800 3rd Qu.: 6.300
## Max. :9.600 Max. :0.6550 Max. :0.6000 Max. :17.500
##
## chlorides free sulfur dioxide total sulfur dioxide
## Min. :0.00900 Min. :0.00100 Min. :0.00600
## 1st Qu.:0.04100 1st Qu.:0.01300 1st Qu.:0.04300
## Median :0.05100 Median :0.02400 Median :0.09900
## Mean :0.05629 Mean :0.02647 Mean :0.09844
## 3rd Qu.:0.07400 3rd Qu.:0.03661 3rd Qu.:0.14300
## Max. :0.10500 Max. :0.07700 Max. :0.27200
##
## density pH sulphates alcohol
## Min. :0.9871 Min. :2.800 Min. :0.2300 Min. : 8.00
## 1st Qu.:0.9931 1st Qu.:3.130 1st Qu.:0.4500 1st Qu.: 9.50
## Median :0.9957 Median :3.230 Median :0.5300 Median :10.22
## Mean :0.9952 Mean :3.235 Mean :0.5369 Mean :10.45
## 3rd Qu.:0.9973 3rd Qu.:3.340 3rd Qu.:0.6100 3rd Qu.:11.20
## Max. :1.0037 Max. :3.630 Max. :0.8500 Max. :14.00
##
## quality color quality class
## 3: 56 red :4797 low :4574
## 4: 373 white:6396 high:6619
## 5:4136
## 6:4760
## 7:1606
## 8: 256
## 9: 6
Y la desviación típica de las variables continuas:
apply(vinos[,c(1:11)], 2, sd)
## fixed acidity volatile acidity citric acid
## 0.889339554 0.129938907 0.113008725
## residual sugar chlorides free sulfur dioxide
## 3.964534828 0.019751699 0.016086736
## total sulfur dioxide density pH
## 0.058921160 0.002844183 0.150358686
## sulphates alcohol
## 0.116513522 1.121871257
Visualizamos las variables (en histograma / gráfico de barras):
col <- c("fixed acidity", "volatile acidity", "citric acid", "residual sugar", "chlorides", "free sulfur dioxide", "total sulfur dioxide" , "density", "pH", "sulphates", "alcohol")
par(mfrow=c(2,3))
for (name in col) {
hist(vinos[,name], prob=TRUE, xlab=name, main = name)
lines(density(na.omit(vinos[,name])))
}
barplot(table(vinos$quality), xlab="quality", main = "quality")
Se observa que la distribucion de la mayoría de las variables es bastante sesgada y parece visualmente a la F-distribution. Suponemos que eso es debido a la natiraleza de indicadores físico-químicos. Ademas, la variable residual sugar tiene una cola por la derecha siquiera despues de tratamiento de los outliers, y su desviación estándar también es relativamente alta.
Hemos visto las distribuciones univariantes del dataset, ahora nos interesa ver las relaciones entre variables, correlaciones, cuáles de las características son más representativas en cuanto a la calidad / color del vino, qué explica la variebilidad natural del dataset y/u otras inferencias que podemos obtener de visualizaciones.
Veamos primero la distribución y las frecuencias relativas de variables por color:
p1 <- ggplot(data=vinos,aes(x=`fixed acidity`,fill=color))+geom_histogram()
p2 <- ggplot(data=vinos,aes(x=`volatile acidity`,fill=color))+geom_histogram()
p3 <- ggplot(data=vinos,aes(x=`citric acid`,fill=color))+geom_histogram()
p4 <- ggplot(data=vinos,aes(x=`residual sugar`,fill=color))+geom_histogram()
p5 <- ggplot(data=vinos,aes(x=`chlorides`,fill=color))+geom_histogram()
p6 <- ggplot(data=vinos,aes(x=`free sulfur dioxide`,fill=color))+geom_histogram()
p7 <- ggplot(data=vinos,aes(x=`total sulfur dioxide`,fill=color))+geom_histogram()
p8 <- ggplot(data=vinos,aes(x=`density`,fill=color))+geom_histogram()+theme(axis.text.x = element_text(angle = 45, hjust = 1))
p9 <- ggplot(data=vinos,aes(x=`pH`,fill=color))+geom_histogram()
p10 <- ggplot(data=vinos,aes(x=`sulphates`,fill=color))+geom_histogram()
p11 <- ggplot(data=vinos,aes(x=`alcohol`,fill=color))+geom_histogram()
p12 <- ggplot(data=vinos,aes(x=`quality`,fill=color))+geom_bar()
p1 + p2 + p3 + p4 + p5 + p6 + p7 + p8 + p9 + p10 + p11 + p12 + plot_layout(ncol = 3) + plot_layout(guides = 'collect')
Se ve que en la mayoría de variables la distribución de vinos tintos y vinos blancos es relativamente diferente. En ‘free sulfur dioxide’ y ‘total sulfur dioxide’ la distribución de tintos tiende bastante más a la izquierda (valores menores) y en otras características hay sesgo a la derecha (valores más altos) - densidad, sulfatos, acidez volatil/fija, cloruros. Y algunas no tienen un sesgo muy pronunciado (alcohol, residual sugar). Muchas caractesristicas bien deben de depender del color del vino.
Como bien tenemos dos conjuntos distintos por el color de vino, podemos comparar las cualidades de los dos, por ejemplo si existe una relacion entre el color de vino y su calidad en el dataset (o bien si son independientes) o por otro lado si hay diferencia en alcohol para vinos tintos/blancos y vinos de calidad baja/alta -con lo que podemos obtener conclusiones tangibles para explicar la naturaleza de vinos.
Hacemos los subconjuntos de datos por el color y por la calidad:
red <- vinos[vinos$color=='red',]
white <- vinos[vinos$color=='white',]
high <- vinos[vinos$`quality class`=='high',]
low <- vinos[vinos$`quality class`=='low',]
Como calidad es la variable discreta, realizamos el test chi-squared para poder hacer inferencia sobre la calidad de los vinos, con la hipótesis nula que la calidad y color son estadísticamente independientes y la hipotesis alternativa que existe una relación entre el color y la calidad.
# Tabla de contingencia
tc<-table( vinos$`quality class`, vinos$color )
tc
##
## red white
## low 2367 2207
## high 2430 4189
# Test
chisq.test(tc, correct=FALSE)
##
## Pearson's Chi-squared test
##
## data: tc
## X-squared = 249.72, df = 1, p-value < 2.2e-16
Segun el p-value podemos concluir que sí hay relacion entre el color de vino y la calidad en los datos observados, pudiendo concluir que los vinos blancos generalmente tienen mejor calidad.
Para comprobar los resultados anteriores, planteamos el contraste para determinar si la proporción de vinos de alta calidad es igual para vinos tintos y blancos. Puesto que tenemos una muestra lo suficinetemente grande, podemos también realizar un test parametrico sobre la proporción.
pR <- dim(vinos[which(vinos$color=='red' & vinos$`quality class`=='high'),])[1]/dim(vinos[vinos$color=='red',])[1];
pW <- dim(vinos[which(vinos$color=='white' & vinos$`quality class`=='high'),])[1]/dim(vinos[vinos$color=='white',])[1];
nR <- dim(vinos[vinos$color=='red',])[1]
nW <- dim(vinos[vinos$color=='white',])[1]
success<-c(pR*nR, pW*nW) # vector de casos de "exito"
nn<-c(nR,nW) # vector de tamaño de muestras
prop.test(success, nn, alternative="two.sided", conf.level=0.95, correct=FALSE)
##
## 2-sample test for equality of proportions without continuity
## correction
##
## data: success out of nn
## X-squared = 249.72, df = 1, p-value < 2.2e-16
## alternative hypothesis: two.sided
## 95 percent confidence interval:
## -0.1667015 -0.1300464
## sample estimates:
## prop 1 prop 2
## 0.5065666 0.6549406
Por el p-value no podemos aceptar la hipotesis nula de igualdad de proporciones, por ello concluimos que las proporciones de vinos de alta calidad son distintos para vinos de distinto color, siendo igualmente el color blanco el que tiende a tener calidad más alta.
Test de hipotesis de contraste sobre la igualdad estadística de la media de alcohol en dos muestras respectivas es un test parametrico si podemos comprobar la suposición de normalidad de variable y segun la igualidad de varianzas eligiremos el test a realizar.
Test de normalidad Shapiro-Wilk (que es mas robusto que el test de Kolmogorov-Smirnov)
shapiro.test(red$alcohol)
##
## Shapiro-Wilk normality test
##
## data: red$alcohol
## W = 0.9347, p-value < 2.2e-16
shapiro.test(sample(white$alcohol,5000))
##
## Shapiro-Wilk normality test
##
## data: sample(white$alcohol, 5000)
## W = 0.95632, p-value < 2.2e-16
No podemos aceptar la normalidad, por ello hacemos la transformación de Box-Cox:
red$alcohol <- BoxCox(red$alcohol, lambda = BoxCoxLambda(red$alcohol))
white$alcohol <- BoxCox(white$alcohol, lambda = BoxCoxLambda(white$alcohol))
shapiro.test(red$alcohol)
##
## Shapiro-Wilk normality test
##
## data: red$alcohol
## W = 0.96551, p-value < 2.2e-16
shapiro.test(sample(white$alcohol,5000))
##
## Shapiro-Wilk normality test
##
## data: sample(white$alcohol, 5000)
## W = 0.97241, p-value < 2.2e-16
A pesar de realizacion de transformacion de Box-Cox, rechazamos la hipotesis nula de normalidad de distribucion. Veamos las gráficas Q-Q:
par(mfrow=c(1,2))
qqnorm(red$alcohol, main='red')
qqline(red$alcohol)
qqnorm(white$alcohol, main='white')
qqline(white$alcohol)
Aunque rechazamos la hipótesis nula de normalidad de datos, con el tamaño de la muestra por el teorema del límite central podemos asumir la normalidad de distribución de alcohol. Con ello realizamos el test de igualdad de varianzas de las muestras:
var.test(red$alcohol, white$alcohol, alternative = "two.sided")
##
## F test to compare two variances
##
## data: red$alcohol and white$alcohol
## F = 0.5761, num df = 4796, denom df = 6395, p-value < 2.2e-16
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
## 0.5464389 0.6074846
## sample estimates:
## ratio of variances
## 0.5760953
De la misma manera, rechazamos la hipótesis nula de igualdad de varianzas por lo que el test estadístico que se realiza es de dos muestras con varianza poblacional desconocida pero distinta.
t.test(red$alcohol, white$alcohol, alternative='two.sided',var.equal=FALSE)
##
## Welch Two Sample t-test
##
## data: red$alcohol and white$alcohol
## t = -4.6589, df = 11189, p-value = 3.216e-06
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -0.0012054766 -0.0004914911
## sample estimates:
## mean of x mean of y
## 0.9028076 0.9036560
Por los resultados del test, con la confianza de 95% aceptamos la hipótesis de que el nivel promedio del alcohol es estadísticamente distinto en vinos blancos y tintos, con la diferencia real entre las dos medias en intervalo [-0.0015334921,-0.0008158837].
Seguimos los mismos pasos para los subconjuntos según la calidad:
Test de normalidad Shapiro-Wilk:
shapiro.test(sample(high$alcohol,5000))
##
## Shapiro-Wilk normality test
##
## data: sample(high$alcohol, 5000)
## W = 0.97549, p-value < 2.2e-16
shapiro.test(low$alcohol)
##
## Shapiro-Wilk normality test
##
## data: low$alcohol
## W = 0.93881, p-value < 2.2e-16
Transformación de Box-Cox:
high$alcohol <- BoxCox(high$alcohol, lambda = BoxCoxLambda(high$alcohol))
low$alcohol <- BoxCox(low$alcohol, lambda = BoxCoxLambda(low$alcohol))
shapiro.test(sample(high$alcohol,5000))
##
## Shapiro-Wilk normality test
##
## data: sample(high$alcohol, 5000)
## W = 0.98265, p-value < 2.2e-16
shapiro.test(low$alcohol)
##
## Shapiro-Wilk normality test
##
## data: low$alcohol
## W = 0.97424, p-value < 2.2e-16
A pesar de realizacion de transformación de Box-Cox, rechazamos la hipótesis nula de normalidad de distribución. Veamos las gráficas Q-Q:
par(mfrow=c(1,2))
qqnorm(high$alcohol, main='high')
qqline(high$alcohol)
qqnorm(low$alcohol, main='low')
qqline(low$alcohol)
Resultado es parecido que en los subconjuntos anteriores, por ello de la misma manera aplicamos el teorema del límite central para asumir la normalidad.
var.test(high$alcohol, low$alcohol, alternative = "two.sided")
##
## F test to compare two variances
##
## data: high$alcohol and low$alcohol
## F = 1.4737, num df = 6618, denom df = 4573, p-value < 2.2e-16
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
## 1.397025 1.554196
## sample estimates:
## ratio of variances
## 1.47371
De la misma manera, rechazamos la hipotesis nula de igualdad de varianzas por lo que el test estadístico que se realiza es de dos muestras con varianza poblacional desconocida pero distinta.
t.test(high$alcohol, low$alcohol, alternative='two.sided',var.equal=FALSE)
##
## Welch Two Sample t-test
##
## data: high$alcohol and low$alcohol
## t = 36.885, df = 10853, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 0.006094820 0.006778973
## sample estimates:
## mean of x mean of y
## 0.9059228 0.8994859
Por los resultados del test, igualmente, con la confianza de 95% aceptamos la hipótesis de que el valor promedio del alcohol es estadísticamente distinto en vinos de calidad alta y baja, con los vinos de alta calidad teniendo el valor promedio de alcohol mas alto.
Por ello, podemos concluir que la variable alcohol puede ser importante tanto para la calidad como para el color del vino, y que hay además una relación estadística entre la calidad y el color.
Podemos ver el correlograma visualizando las correlaciones entre las variables:
d <- vinos[1:11]
d$quality = as.numeric(vinos$quality)
d$color = as.numeric(vinos$color)
M<-cor(d)
corrplot(M, method="ellipse", type='lower',tl.col="black", tl.srt=45)
Se observa que, muchos atributos tienen correlación entre sí. Por ejemplo, fuerte correlación positiva presentan total sulfur dioxide y free sulfur dioxide, esperadamente. En general, density tiene correlación positiva con casi todas las características (sin contar color y calidad), bastante alta con “residual sugar”, y solo con “alcohol” está fuerte y negativamente correlacionada. Las correlaciones entre atributos químicos pueden ser explicados por adición de las sustancias con el fin de nivelar/acentuar las otras características .
Por otro lado, calidad no parece tener correlaciones muy importantes positivas o negativas con muchas de las características, no obstante las que mâs le puedan influir son alcohol (positivamente) y densidad (negativamente).
No obstante, el color está más correlacionado con los atributos físico-químicos, como con “total sulfur dioxide” (positivamente), “volatile acidity” (negativamente), etc, donde los signos de correlación representaran: negativo -> vino tinto, positivo -> vino blanco. Hay que tener en cuenta que la dimensionalidad de la parte de vinos tintos en el dataset es distinta a la de vinos blancos, aunque debe de ser suficiente para obtener las correlaciones correctas.
Algunas de correlaciones más signicativas:
p1 <- ggplot(data = vinos,aes(x=`free sulfur dioxide`,y=`total sulfur dioxide`))+geom_jitter(color="#0072B2", alpha=0.3) + geom_smooth(method='lm')
p2 <- ggplot(data = vinos,aes(x=`density`,y=`residual sugar`))+geom_jitter(color="#0072B2", alpha=0.3) + geom_smooth(method='lm')
p3 <- ggplot(data = vinos,aes(x=`citric acid`,y=`volatile acidity`)) + geom_jitter(color="#FF9966", alpha=0.3) + geom_smooth(method='lm')
p4 <- ggplot(data = vinos,aes(x=`alcohol`,y=`residual sugar`)) + geom_jitter(color="#FF9966", alpha=0.3) + geom_smooth(method='lm')
p1 + p2 + p3 + p4
Es inetersante ver qué color tiene las correlaciones más fuertes con atributos físicos y químicos de vinos que calidad.
Veamos algunas de las características:
p1 <- ggplot(data = vinos,aes(x=`volatile acidity`,fill=color))+geom_density(alpha=0.7, outline.type = 'lower')
p2 <- ggplot(data = vinos,aes(x=`total sulfur dioxide`,fill=color))+geom_density(alpha=0.7, outline.type = 'lower')
p3 <- ggplot(data = vinos,aes(x=`chlorides`,fill=color))+geom_density(alpha=0.7, outline.type = 'lower')
p4 <- ggplot(data = vinos,aes(x=`sulphates`,fill=color))+geom_density(alpha=0.7, outline.type = 'lower')
p1 + p2 + p3 + p4 + plot_layout(guides = 'collect')
Las dos primeras características (“total sulphur dioxide” y “volatile acidity”) resultan en distribuciones de densidad muy distintas entre tipos de vinos que confirma la correlación.
Además, podemos ver las agrupaciones por color visualizando dos variables opuestamente correlacionadas:
ggplot(data = vinos, aes(x = `total sulfur dioxide`, y = `volatile acidity`, color = `color`)) + geom_point(alpha = 0.3)
Hay cierta tendencia bastante definida visualmente de agrupaciones de vinos tintos/blancos aunque las observaciones están todavía entremezcladas.
Fijémonos en estas dos características:
p1 <- ggplot(data = vinos,
aes(x = `total sulfur dioxide`, y = `volatile acidity`, color = `color`)) +
geom_point(alpha = 0.3)+facet_wrap(~color)
p2 <- ggplot(data = vinos,
aes(x = `total sulfur dioxide`, y = `volatile acidity`, color = `color`)) +
geom_point(alpha = 0.3)+facet_wrap(~`quality class`)
p1 + p2 + plot_layout(guides = 'collect') + plot_layout(ncol = 1)
Las agrupaciones se mantienen en distintos grupos de calidad de vino.
# Total sulfur dioxide
p1 <- ggplot(data = vinos,aes(x=`total sulfur dioxide`,fill=`color`))+geom_density(alpha=0.5, outline.type = 'lower')+facet_wrap(~color)
p2 <- ggplot(data = vinos,aes(x=`quality class`, y=`total sulfur dioxide`, color = `color`))+geom_boxplot()+facet_wrap(~`color`)
p1 + p2 + plot_layout(guides = 'collect')
# Volatile acidity
p1 <- ggplot(data = vinos,aes(x=`volatile acidity`,fill=`color`))+geom_density(alpha=0.5, outline.type = 'lower')+facet_wrap(~color)
p2 <- ggplot(data = vinos,aes(x=`quality class`, y=`volatile acidity`, color = `color`))+geom_boxplot()+facet_wrap(~`color`)
p1 + p2 + plot_layout(guides = 'collect')
Se ven claramente las diferencias en distribuciones de ambas características por color de vino y correlaciones positiva y negativa. Ya que “volatile acidity” también presntaba correlación con “quality”, se observan diferencias por cada clase de calidad, pero principalmente en vinos tintos.
La densidad y el alcohol pueden ser las características del vino que son las más responsables por la calidad del vino del dataset. También, hay correlaciones más bajas con volatile acidity, sin embargo, como hemos visto, principalmente se manifestan en un tipo de vinos (vinos tintos).
Veamos las dos variables con más detalle para ver la relación entre ellas:
ggplot(data = vinos, aes(x = density, y = alcohol, color = `quality class`)) + geom_point(alpha = 0.3)
p1 <- ggplot(data = vinos,
aes(x = density, y = alcohol, color = `quality class`)) +
geom_point(alpha = 0.3)+facet_wrap(~color)
p2 <- ggplot(data = vinos,
aes(x = density, y = alcohol, color = `quality class`)) +
geom_point(alpha = 0.3)+facet_wrap(~`quality class`)
p1 + p2 + plot_layout(guides = 'collect') + plot_layout(ncol = 1)
Aunque las clases están entremezclados, también se ve la tendencia de agrupaciones.
# Alcohol
p1 <- ggplot(data = vinos,aes(x=alcohol,fill=`quality class`))+geom_density(alpha=0.5, outline.type = 'lower')+facet_wrap(~`quality class`)
p2 <- ggplot(data = vinos,aes(x=color, y=alcohol, color = `quality class`))+geom_boxplot()+facet_wrap(~`quality class`)
p1 + p2 + plot_layout(guides = 'collect')
# Density
p1 <- ggplot(data = vinos,aes(x=density,fill=`quality class`))+geom_density(alpha=0.5, outline.type = 'lower')+facet_wrap(~`quality class`)+theme(axis.text.x = element_text(angle = 45, hjust = 1))
p2 <- ggplot(data = vinos,aes(x=color, y=density, color = `quality class`))+geom_boxplot()+facet_wrap(~`quality class`)
p1 + p2 + plot_layout(guides = 'collect')
Observamos las correlaciones (positiva con alcohol y negativa con densidad) en ambos tipos/colores de vinos, aunque la densidad es la característica con una influencia más fuerte en vinos blancos. Se visualiza bien cómo el nivel de alcohol y densidad cambian en distintas particiones de la calidad del vino del dataset.
Hay tendencias de agrupamiento por valores de nivel de alcohol y densidad, pero al final la calidad debe de depender de un conjunto más grande de características físico-químicas ya que dos variables no explican totalmente la calidad, ni calidad explica la variebilidad del dataset.
Finalmente, para aproximarnos al conjunto de características, podemos visualizar la relación entre 4 variables: density, alcohol, chlorides y quality class. Chlorides es la siguiente característica que tenía correlacion fuerte con quality despues de volatile acidity que lo presentaba mayormente en vinos tintos.
ggplot(data = vinos, aes(x = density, y = alcohol, color = `quality class`)) + geom_point(aes(size = `chlorides`), alpha = 0.3)
No hay relación lineal, pero se explican algunas observaciones “alejadas” (pertenecen a una partición de calidad pero se visualizan donde prevalece otra partición, pues tienen chlorides, por ejemplo, bajos, frente a la mayoría de su partición. Al final y a cabo muchos de los atributos influen a la densidad).
Como hemos observado, no hay una relacion lineal simple entre las variables, por ello se puede intentar explicar la calidad y el color linealmente con varios regresores a traves de algoritmo de regresion logística, dado que tanto el color como la clase de calidad son variables dicotómicas.
model1 <- glm(`quality class`~`fixed acidity`+`volatile acidity`+`citric acid`+`residual sugar`+chlorides+`free sulfur dioxide`+`total sulfur dioxide`+density+pH+sulphates+alcohol, vinos, family=binomial(link=logit))
summary(model1)
##
## Call:
## glm(formula = `quality class` ~ `fixed acidity` + `volatile acidity` +
## `citric acid` + `residual sugar` + chlorides + `free sulfur dioxide` +
## `total sulfur dioxide` + density + pH + sulphates + alcohol,
## family = binomial(link = logit), data = vinos)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.4331 -1.0782 0.5563 0.9612 1.8516
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.051e+02 1.602e+01 6.561 5.34e-11 ***
## `fixed acidity` 3.658e-02 2.728e-02 1.341 0.180
## `volatile acidity` -2.224e+00 1.924e-01 -11.561 < 2e-16 ***
## `citric acid` 1.710e-01 2.023e-01 0.845 0.398
## `residual sugar` 7.823e-02 8.227e-03 9.509 < 2e-16 ***
## chlorides -2.164e+00 1.657e+00 -1.306 0.192
## `free sulfur dioxide` 1.510e+01 2.113e+00 7.145 9.01e-13 ***
## `total sulfur dioxide` -4.489e+00 6.633e-01 -6.769 1.30e-11 ***
## density -1.115e+02 1.609e+01 -6.928 4.27e-12 ***
## pH -6.574e-02 1.653e-01 -0.398 0.691
## sulphates 1.925e+00 2.241e-01 8.591 < 2e-16 ***
## alcohol 5.487e-01 3.098e-02 17.713 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 15141 on 11192 degrees of freedom
## Residual deviance: 13308 on 11181 degrees of freedom
## AIC: 13332
##
## Number of Fisher Scoring iterations: 4
Según el modelo, las variables fixed acidity, citric acid, pH no tienen significancia estadística para el modelo (p.value > 0.05 por lo que aceptamos la hipótesis que significancia es igual a cero).
Para ver los odds ratio para cada unidad de las caracteristicas:
exp(coefficients(model1))
## (Intercept) `fixed acidity` `volatile acidity`
## 4.340860e+45 1.037258e+00 1.081545e-01
## `citric acid` `residual sugar` chlorides
## 1.186468e+00 1.081376e+00 1.148246e-01
## `free sulfur dioxide` `total sulfur dioxide` density
## 3.614470e+06 1.122705e-02 3.898240e-49
## pH sulphates alcohol
## 9.363707e-01 6.857162e+00 1.730965e+00
Puesto que las caracteristicas con valores superiores a uno son alcohol y residual sugar, podemos considerar alcohol el más influyente a la probabilidad en el modelo obtenido.
La bondad de ajuste se obtiene a traves del índice de Akaike AIC, que en este caso asciende a 13336, un valor muy alto. Podremos compararlo con el modelo posterior de regresión sobre el color.
model2 <- glm(color~`fixed acidity`+`volatile acidity`+`citric acid`+`residual sugar`+chlorides+`free sulfur dioxide`+`total sulfur dioxide`+density+pH+sulphates+alcohol, vinos, family=binomial(link=logit))
summary(model2)
##
## Call:
## glm(formula = color ~ `fixed acidity` + `volatile acidity` +
## `citric acid` + `residual sugar` + chlorides + `free sulfur dioxide` +
## `total sulfur dioxide` + density + pH + sulphates + alcohol,
## family = binomial(link = logit), data = vinos)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -4.0193 -0.0167 0.0062 0.0407 3.0691
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.023e+03 5.770e+01 17.722 < 2e-16 ***
## `fixed acidity` -6.105e-01 1.148e-01 -5.317 1.05e-07 ***
## `volatile acidity` -8.871e+00 8.202e-01 -10.816 < 2e-16 ***
## `citric acid` 3.464e+00 9.788e-01 3.539 0.000401 ***
## `residual sugar` 4.717e-01 4.297e-02 10.979 < 2e-16 ***
## chlorides -8.449e+01 6.535e+00 -12.930 < 2e-16 ***
## `free sulfur dioxide` -7.656e+01 1.055e+01 -7.255 4.02e-13 ***
## `total sulfur dioxide` 7.722e+01 4.060e+00 19.017 < 2e-16 ***
## density -1.001e+03 5.749e+01 -17.421 < 2e-16 ***
## pH -1.466e+00 6.957e-01 -2.107 0.035087 *
## sulphates -8.904e+00 9.658e-01 -9.220 < 2e-16 ***
## alcohol -1.162e+00 1.270e-01 -9.147 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 15287.58 on 11192 degrees of freedom
## Residual deviance: 898.12 on 11181 degrees of freedom
## AIC: 922.12
##
## Number of Fisher Scoring iterations: 9
Se observan en este caso todas las variables significativas para el modelo, con los odds ratio respectivos de:
exp(coefficients(model2))
## (Intercept) `fixed acidity` `volatile acidity`
## Inf 5.430872e-01 1.404404e-04
## `citric acid` `residual sugar` chlorides
## 3.194939e+01 1.602773e+00 2.023956e-37
## `free sulfur dioxide` `total sulfur dioxide` density
## 5.627164e-34 3.429002e+33 0.000000e+00
## pH sulphates alcohol
## 2.308147e-01 1.358281e-04 3.130077e-01
según los que las variables explicativas parecen ser total sulfur dioxide, alcohol.
El modelo es mas explicativo que el anterior, puesto que la bondad de ajuste es mejor con el AIC de 823.
Podemos decir que no es facil encontrar un modelo lineal multiple que explique la varianza del dataset tanto para la probabilidad del color como de la calidad, entonces que el dataset no es linealmente separable. Otra manera poder explicar la varianza es realizar un estudio no supervisado y analizar las agrupaciones resultantes.
Normalizamos el dataset (solo los features), puesto que el algoritmo de clustering se basa en las distancias entre las observaciones.
vinos_norm <- scale(vinos[,1:11])
Tal y como hemos dicho anteriormente, aunque hemos estumado el numero de intervalos de calidad de 2 (alta/baja), no tiene por que ser un número óptimo de clusters segun las caracteristicas físico-químicas, por lo que estudiaremos las agrupaciones sin referencia a la marca de calidad y las aproximaremos al dataset con atributos caractrísticos del vino.
Con la librería factoextra y la función fviz_nbclust podemos visualizar el número óptimo calculado con metodos por siluetas medias y wss (suma de cuadrados):
# Metodo por siluetas medias
fviz_nbclust(vinos_norm, kmeans, method = "silhouette")
# Metodo wss
fviz_nbclust(vinos_norm, kmeans, method = "wss")
El número óptimo de clusters parece ser 2, lo que está más claro por el metodo de siluetas pero no por el metodo de “codo”, como se ha preliminarmente definido para el número de particiones de calidad.
set.seed(1234)
fit <- kmeans(vinos_norm, 2)
# Tamaño de los clusters obtenidos
fit$size
## [1] 4807 6386
Visualizamos los clusters:
fviz_cluster(fit, geom = "point", data = vinos_norm)
Se observan que están bastante poco entremezclados, siendo el segundo cluster más poblado.
Se puede ver la feature importance de las dimensiones de la gráfica que serán las mismas que los 2 primeros componentes de PCA:
var <- get_pca_var(prcomp(vinos_norm))
p1 <- fviz_contrib(prcomp(vinos_norm), choice = "var", axes = 1, top = 5)
p2 <- fviz_contrib(prcomp(vinos_norm), choice = "var", axes = 2, top = 5)
p1 + p2
Para poder aproximar las agrupaciones al dataset, visualizamos un gráfico con los atributos explicativos del dataset (no componentes).
Para visualizar los gráficos con los atributos total sulfur dioxide y chlorides (dimensión 1):
ggplot(vinos,aes(`total sulfur dioxide`, `chlorides`, color=as.factor(fit$cluster)))+geom_point()
ggplot(vinos,aes(`total sulfur dioxide`, `chlorides`, color=`quality class`))+geom_point()
ggplot(vinos,aes(`total sulfur dioxide`, `chlorides`, color=color))+geom_point()
Para visualizar ahora los gráficos con los atributos density y alcohol (dimensión 2):
ggplot(vinos,aes(`density`, `alcohol`, color=as.factor(fit$cluster)))+geom_point()
ggplot(vinos,aes(`density`, `alcohol`, `chlorides`, color=`quality class`))+geom_point()
ggplot(vinos,aes(`density`, `alcohol`, `chlorides`, color=color))+geom_point()
En este caso la aproximación de clusters puede captar la natiraleza de las variables segun su color pero no segun quality class, lo que significa que la discretización realizada no es muy significativa, pero que las caracteristicas tienen que ser distintas segun el color.
A partir de los resultados obtenidos, ¿cuáles son las conclusiones? ¿Los resultados permiten responder al problema?
En principio no podemos concluir que sea posible determinar una calidad de vino segun su sensory data con alta precisión. Tenemos evidencia que el color sí que esta correlacionado con los atributos.
(no acabado)